Tree Density from LiDAR

Author

Jay Matsushiba

Published

April 25, 2025

Setup

Code
library(bowen.biodiversity.webapp)
library(here)
library(sf)
library(leaflet)
library(leaflet.extras)
library(future)
library(dplyr)

# install.packages('lasR', repos = 'https://r-lidar.r-universe.dev')
library(lasR)

Extract Tree Tops from LiDAR Data

LiDAR data files were downloaded from LidarBC

The lasR package was used to extract tree tops data from the LiDAR data.

Code
set_parallel_strategy(concurrent_files(4))
# Read las 
ctg <- readLAScatalog(here("temp/las"))

del = triangulate(filter = keep_ground())
norm = transform_with(del, "-")
write = write_las(here("temp/normalized/*_normalized.las"))
chm = rasterize(2, "max", ofile = here("temp/chm/chm.tif"))
lm = local_maximum_raster(chm, 5)
pipeline2 = reader() + del + norm + write + chm + lm

ans = exec(pipeline2, on = ctg)
ttops <- ans$local_maximum %>% st_transform(4326)
st_write(ttops, here("data-raw/ttops.gpkg"))

Visualize Tree Tops Data with Heatmap

Code
ttops <- st_read(here("data-raw/ttops.gpkg")) 
Reading layer `ttops' from data source 
  `/home/jay/Programming_Projects/bowen.biodiversity.webapp/data-raw/ttops.gpkg' 
  using driver `GPKG'
Simple feature collection with 660834 features and 0 fields
Geometry type: POINT
Dimension:     XYZ
Bounding box:  xmin: -123.4483 ymin: 49.32728 xmax: -123.3 ymax: 49.42478
z_range:       zmin: 2 zmax: 1383.31
Geodetic CRS:  WGS 84
Code
ttops_coord <- ttops %>% st_coordinates() %>% as.data.frame()
leaflet(ttops_coord) %>%      
  addProviderTiles(providers$CartoDB.Positron) %>%
  addHeatmap(lng=~X, lat=~Y, blur= 15, max = 65, radius = 30)